home *** CD-ROM | disk | FTP | other *** search
/ Super PC 34 / Super PC 34 (Shareware).iso / spc / UTIL / DJGPP2 / V2 / DJTST200.ZIP / tests / libc / ansi / math / elefunt / texp.c < prev    next >
Encoding:
C/C++ Source or Header  |  1995-06-01  |  6.1 KB  |  262 lines
  1. /* -*-C-*- texp.c */
  2.  
  3. #include "elefunt.h"
  4.  
  5. /*
  6. float EXP(float);
  7.  
  8. float EXP(float x)
  9. {
  10.     if (x < 0.0)
  11.         return(1.0/exp(-x));
  12.     return(exp(x));
  13. }
  14. */
  15.  
  16. #define EXP exp
  17.  
  18. /***********************************************************************
  19. #     program to test exp
  20. #
  21. #     data required
  22. #
  23. #        none
  24. #
  25. #     subprograms required from this package
  26. #
  27. #        machar - an environmental inquiry program providing
  28. #                 information on the floating-point arithmetic
  29. #                 system.  note that the call to machar can
  30. #                 be deleted provided the following four
  31. #                 parameters are assigned the values indicated
  32. #
  33. #                 ibeta - the radix of the floating-point system
  34. #                 it    - the number of base-ibeta digits in the
  35. #                         significand of a floating-point number
  36. #                 xmin  - the smallest non-vanishing floating-point
  37. #                         power of the radix
  38. #                 xmax  - the largest finite floating-point no.
  39. #
  40. #        ran(k) - a function subprogram returning random real
  41. #                 numbers uniformly distributed over (0,1)
  42. #
  43. #
  44. #     standard fortran subprograms required
  45. #
  46. #         abs, aint, alog, AMAX1, exp, float, sqrt
  47. #
  48. #
  49. #     latest revision - december 6, 1979
  50. #
  51. #     author - w. j. cody
  52. #              argonne national laboratory
  53. #
  54. #
  55. ***********************************************************************/
  56.  
  57. void
  58. texp()
  59. {
  60.     int i,
  61.         ibeta,
  62.         iexp,
  63.         irnd,
  64.         it,
  65.         i1,
  66.         j,
  67.         k1,
  68.         k2,
  69.         k3,
  70.     machep,
  71.         maxexp,
  72.         minexp,
  73.         n,
  74.         negep,
  75.         ngrd;
  76.     float a,
  77.         ait,
  78.         albeta,
  79.         b,
  80.         beta,
  81.         d,
  82.         del,
  83.         eps,
  84.         epsneg,
  85.         r6,
  86.         r7,
  87.         v,
  88.         w,
  89.     x,
  90.         xl,
  91.         xmax,
  92.         xmin,
  93.         xn,
  94.         x1,
  95.         y,
  96.         z,
  97.         zz;
  98.  
  99.     machar(&ibeta, &it, &irnd, &ngrd, &machep, &negep, &iexp, &minexp,
  100.         &maxexp, &eps, &epsneg, &xmin, &xmax);
  101.     beta = (float) ibeta;
  102.     albeta = ALOG(beta);
  103.     ait = (float) it;
  104.     v = 0.0625e0L;
  105.     a = TWO;
  106.     b = ALOG(a) * 0.5e0L;
  107.     a = -b + v;
  108.     d = ALOG(0.9e0L * xmax);
  109.     n = 2000;
  110.     xn = (float) n;
  111.     i1 = 0;
  112.  
  113.     /* random argument accuracy tests */
  114.     for (j = 1; j <= 3; ++j)
  115.     {
  116.     k1 = 0;
  117.     k3 = 0;
  118.     x1 = ZERO;
  119.     r6 = ZERO;
  120.     r7 = ZERO;
  121.     del = (b - a) / xn;
  122.     xl = a;
  123.  
  124.     for (i = 1; i <= n; ++i)
  125.     {
  126.         x = del * ran(i1) + xl;
  127.  
  128.         /* purify arguments */
  129.         y = x - v;
  130.         if (y < ZERO)
  131.         x = y + v;
  132.         z = EXP(x);
  133.         zz = EXP(y);
  134.         if (j == 1)
  135.         z = z - z * 6.058693718652421388e-2L;
  136.         else
  137.         {
  138.         if (ibeta != 10)
  139.             /* z = z *.0625e0L - z * 2.4453321046920570389e-3L; */
  140.             z = z * EXP(-2.8125L);
  141.         else
  142.             z = z * 6.0e-2L + z * 5.466789530794296106e-5L;
  143.         }
  144.         w = ONE;
  145.         if (zz != ZERO)
  146.         w = (z - zz) / zz;
  147.         if (w < ZERO)
  148.         k1 = k1 + 1;
  149.         if (w > ZERO)
  150.         k3 = k3 + 1;
  151.         w = ABS(w);
  152.         if (w > r6)
  153.         {
  154.         r6 = w;
  155.         x1 = x;
  156.         }
  157.         r7 = r7 + w * w;
  158.         xl = xl + del;
  159.     }
  160.  
  161.     k2 = n - k3 - k1;
  162.     r7 = sqrt(r7 / xn);
  163.  
  164.     printf("\fTEST OF EXP(X-" F7P4F ") VS EXP(X)/EXP(" F7P4F ")\n\n\n", v, v);
  165.     printf("%7d RANDOM ARGUMENTS WERE TESTED FROM THE INTERVAL\n", n);
  166.     printf("      (" F15P4E "," F15P4E ")\n\n\n", a, b);
  167.     printf(" EXP(X-V) WAS LARGER%6d TIMES,\n", k1);
  168.     printf("              AGREED%6d TIMES, AND\n", k2);
  169.     printf("         WAS SMALLER%6d TIMES.\n\n\n", k3);
  170.     printf(
  171. " THERE ARE %4d BASE %4d SIGNIFICANT DIGITS IN A FLOATING-POINT NUMBER\n\n\n",
  172.         it, ibeta);
  173.     w = -999.0e0;
  174.     if (r6 != ZERO)
  175.         w = ALOG(ABS(r6)) / albeta;
  176.     printf(" THE MAXIMUM RELATIVE ERROR OF" F15P4E " = %4d **" F7P2F "\n",
  177.         r6, ibeta, w);
  178.     printf("    OCCURRED FOR X =" F17P6E "\n", x1);
  179.     w = AMAX1(ait + w, ZERO);
  180.     printf(
  181.         " THE ESTIMATED LOSS OF BASE%4d SIGNIFICANT DIGITS IS" F7P2F "\n\n\n",
  182.         ibeta, w);
  183.     w = -999.0e0;
  184.     if (r7 != ZERO)
  185.         w = ALOG(ABS(r7)) / albeta;
  186.     printf(" THE ROOT MEAN SQUARE RELATIVE ERROR WAS" F15P4E " = %4d **" F7P2F "\n",
  187.         r7, ibeta, w);
  188.     w = AMAX1(ait + w, ZERO);
  189.     printf(
  190.         " THE ESTIMATED LOSS OF BASE%4d SIGNIFICANT DIGITS IS" F7P2F "\n\n\n",
  191.         ibeta, w);
  192.     if (j == 2)
  193.     {
  194.         a = -TWO * a;
  195.         b = TEN * a;
  196.         if (b < d)
  197.         b = d;
  198.     }
  199.     else
  200.     {
  201.         v = 45.0e0L / 16.0e0L;
  202.         a = -TEN * b;
  203.         b = 4.0e0L * xmin * ipow((float) beta, it);
  204.         b = ALOG(b);
  205.     }
  206.     }
  207.  
  208.     /* special tests */
  209.     printf("\fSPECIAL TESTS\n\n\n");
  210.     printf(" THE IDENTITY EXP(X)*EXP(-X) = 1.0  WILL BE TESTED.\n\n");
  211.     printf("       X        F(X)*F(-X) - 1\n\n");
  212.  
  213.     for (i = 1; i <= 5; ++i)
  214.     {
  215.     x = ran(i1) * beta;
  216.     y = -x;
  217.     z = EXP(x) * EXP(y) - ONE;
  218.     printf(F15P7E F15P7E "\n\n", x, z);
  219.     }
  220.  
  221.     printf("\n\n TEST OF SPECIAL ARGUMENTS\n\n\n");
  222.     x = ZERO;
  223.     y = EXP(x) - ONE;
  224.     printf(" EXP(0.0) - 1.0E0 = " F15P7E "\n\n", y);
  225.     x = AINT(ALOG(xmin));
  226.     y = EXP(x);
  227.     printf(" EXP(" F13P6E ") =" F13P6E "\n\n", x, y);
  228.     x = AINT(ALOG(xmax));
  229.     y = EXP(x);
  230.     printf(" EXP(" F13P6E ") =" F13P6E "\n\n", x, y);
  231.     x = x / TWO;
  232.     v = x / TWO;
  233.     y = EXP(x);
  234.     z = EXP(v);
  235.     z = z * z;
  236.     printf(" IF EXP(" F13P6E ") = " F13P6E " IS NOT ABOUT\n", x, y);
  237.     printf(" EXP(" F13P6E ")**2 =" F13P6E " THERE IS AN ARG RED ERROR\n", v, z);
  238.  
  239.     /* test of error returns */
  240.  
  241.     printf("\fTEST OF ERROR RETURNS\n\n\n");
  242.     x = -ONE / sqrt(xmin);
  243.     printf(" EXP WILL BE CALLED WITH THE ARGUMENT" F15P4E "\n",x);
  244.     printf(" THIS SHOULD TRIGGER AN ERROR MESSAGE\n\n\n");
  245.     fflush(stdout);
  246.     errno = 0;
  247.     y = EXP(x);
  248.     if (errno)
  249.     perror("exp()");
  250.     printf(" EXP RETURNED THE VALUE" F13P6E "\n\n\n\n", y);
  251.     x = -x;
  252.     printf(" EXP WILL BE CALLED WITH THE ARGUMENT" F13P6E "\n", x);
  253.     printf(" THIS SHOULD TRIGGER AN ERROR MESSAGE\n\n\n");
  254.     fflush(stdout);
  255.     errno = 0;
  256.     y = EXP(x);
  257.     if (errno)
  258.     perror("exp()");
  259.     printf(" EXP RETURNED THE VALUE" F13P6E "\n\n\n\n", y);
  260.     printf(" THIS CONCLUDES THE TESTS\n");
  261. }
  262.